
function beta_hat=est_cauchit_rewt(y,x,tau);

%{
This routine runs a cauchit regression of y on x (based on a linear index). 
y is an (n x 1) vector consisting of zeros and ones, while x is a (n x k) 
matrix of explanatory variables. 

The sample is assumed to be a stratified sample with strata given by Y=0 
and Y=1. The uncoditional probability that Y=1 is tau. The sample proportion
of observations with Y=1 cannot be taken as an estimate of tau. 

The output is a (k x 1) column vector of estimated logit coefficients and 
its (k x k) covariance matrix.
%}

dim=size(x);
k=dim(2); %number of variables in x

index_y0=find(y==0);
index_y1=find(y==1);

x0=x(index_y0,:);
x1=x(index_y1,:);

beta_ini=zeros(k,1); 
options=optimset('Display', 'off'); 
beta_ini=fsolve(@(beta) logit_score_rewt(beta, x0, x1, tau),  beta_ini, options); %logit starting value (logit is robust)

beta_hat=fsolve(@(beta) cauchit_score_rewt(beta, x0, x1, tau),  beta_ini, options);

%Note: the command @(beta) logit_score(beta, y, x) designates logit_score
%as a function of the single argument beta for use with the fsolve command. 
%(Formally, it creates an anonymous function.) Thus, beta_ini is the 
%initial value fsolve will use to search for the zeros of the score fn.

%Log-likelihood:
% L=sum( log( tcdf(x*beta_hat,1) ).*(y==1) + log( tcdf(-x*beta_hat, 1)).*(y==0) );
% beta_hat=beta_hat'; %1 x k

 



